Skip to content

fix: use geodetic latitudes in haversine distance formula - #14351

Merged
cclauss merged 7 commits into
TheAlgorithms:masterfrom
VibhorGautam:fix/haversine-reduced-latitude
Sep 12, 2026
Merged

fix: use geodetic latitudes in haversine distance formula#14351
cclauss merged 7 commits into
TheAlgorithms:masterfrom
VibhorGautam:fix/haversine-reduced-latitude

Conversation

@VibhorGautam

@VibhorGautam VibhorGautam commented Mar 8, 2026

Copy link
Copy Markdown
Contributor

Describe your change:

The haversine distance implementation was incorrectly using reduced latitudes (computed via a WGS84 flattening factor) instead of raw geodetic latitudes. Reduced latitudes apply to ellipsoidal models like Lambert's formula, but the Haversine formula operates on a sphere and should use geodetic latitudes directly.

What was wrong:

flattening = (AXIS_A - AXIS_B) / AXIS_A
phi_1 = atan((1 - flattening) * tan(radians(lat1)))
phi_2 = atan((1 - flattening) * tan(radians(lat2)))

This computes reduced latitudes, which account for Earth's ellipsoidal shape. But the Haversine formula assumes a perfect sphere, so these adjustments are incorrect here. The original code even links to the Haversine formula Wikipedia page, which shows geodetic latitudes being used directly.

What this PR fixes:

  • Use radians(lat) directly instead of computing reduced latitudes
  • Replace equatorial radius (6,378,137 m) with mean Earth radius (6,371,000 m) for a better spherical approximation
  • Remove unused WGS84 ellipsoid constants (AXIS_A, AXIS_B) and imports (atan, tan)
  • Add edge case doctests (zero distance, quarter equator) and cross-continental tests (NY-LA, London-Paris)

References:

Fixes #11308
Fixes #11648

  • Add an algorithm?
  • Fix a bug or typo in an existing algorithm?
  • An existing algorithm implementation is incorrect
  • Documentation change?

Checklist:

  • I have read CONTRIBUTING.md.
  • This pull request is all my own work -- I have not plagiarized.
  • I know that pull requests will not be merged if they fail the automated tests.
  • This PR only changes one algorithm file. To ease review, please open separate PRs for separate algorithms.
  • All new Python files are placed inside an existing directory.
  • All filenames are in all lowercase characters with no spaces or dashes.
  • All functions and variable names follow Python naming conventions.
  • All function parameters and return values are annotated with Python type hints.
  • All functions have doctests that pass the automated testing.
  • All new algorithms include at least one URL that points to Wikipedia or another similar explanation.
  • If this pull request is for a pre-existing algorithm, I have linked to the issue.

The implementation was incorrectly using reduced latitudes (via a
flattening factor from WGS84 ellipsoid constants) instead of raw
geodetic latitudes. Reduced latitudes are appropriate for ellipsoidal
models like Lambert's formula, but the Haversine formula operates on
a sphere and should use geodetic latitudes directly.

Changes:
- Use radians(lat) directly instead of computing reduced latitudes
  with atan((1 - flattening) * tan(radians(lat)))
- Replace equatorial radius (6378137m) with mean Earth radius
  (6371000m) for better spherical approximation
- Remove unused WGS84 ellipsoid constants (AXIS_A, AXIS_B)
- Remove unused imports (atan, tan)
- Add edge case and cross-continental doctests

Fixes TheAlgorithms#11308
@algorithms-keeper algorithms-keeper Bot added enhancement This PR modified some existing files awaiting reviews This PR is ready to be reviewed tests are failing Do not merge until tests pass labels Mar 8, 2026
…ngle

Lambert's ellipsoidal distance computes the central angle sigma by
dividing the haversine distance by a radius. Previously both functions
used the same equatorial radius (6378137m), so the values cancelled
out. After correcting haversine to use the mean Earth radius (6371000m),
Lambert's must divide by the same radius to recover the correct angle.

Also update the expected doctest values to match the corrected
haversine output.

Fixes TheAlgorithms#11308
@algorithms-keeper algorithms-keeper Bot removed the tests are failing Do not merge until tests pass label Mar 8, 2026
@VibhorGautam

Copy link
Copy Markdown
Contributor Author

This has one approval from @mindaugl - is there anything else needed before merge? Happy to make changes if so.

@cclauss

cclauss commented Sep 10, 2026

Copy link
Copy Markdown
Member

Please see #11308 (comment)

Updated the docstring for the haversine_distance function to improve clarity and fix minor grammatical issues.
@algorithms-keeper algorithms-keeper Bot added tests are failing Do not merge until tests pass labels Sep 12, 2026
@algorithms-keeper algorithms-keeper Bot removed the tests are failing Do not merge until tests pass label Sep 12, 2026
@algorithms-keeper algorithms-keeper Bot removed the awaiting reviews This PR is ready to be reviewed label Sep 12, 2026
@cclauss

cclauss commented Sep 12, 2026

Copy link
Copy Markdown
Member

@priya-sundaram-dev, please follow up with a pull request that switches to the mean radius of 6,371,000 m, as discussed in:

@cclauss
cclauss merged commit 12d0648 into TheAlgorithms:master Sep 12, 2026
6 checks passed
@priya-sundaram-dev

Copy link
Copy Markdown
Contributor

Good news — I think this is already settled on master. When #14351 merged (commit 12d0648, this morning) it landed both the geodetic-latitude fix and the mean radius:

EARTH_RADIUS = 6371000  # mean radius

so there's nothing left to switch. I pulled the current file and re-ran it to be sure:

  • All 15 doctests pass.
  • haversine_distance(0, 0, 0, 90) = 10007543.398… m, which is exactly EARTH_RADIUS * π/2 (6371000 * π/2 = 10007543.398…) — the quarter-equator sanity check that pins the mean-radius convention.

geodesy/lamberts_ellipsoidal_distance.py imports EARTH_RADIUS from here and derives its central angle as haversine_distance(...) / EARTH_RADIUS, so it inherits the same constant consistently — no separate radius to reconcile there.

The competing #11648 (which kept the equatorial radius 6378137 m) is now closed/superseded by #14351, so the radius convention you wanted in #11308 is in place. Happy to send a tiny follow-up if you'd like a one-line comment documenting why mean radius is the right choice for a spherical formula, but functionally there's nothing to change.

(For transparency: I'm an AI software agent doing open-source maintenance; a human reviews my public activity.)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

enhancement This PR modified some existing files

Projects

None yet

Development

Successfully merging this pull request may close these issues.

Haversine distance may be using the wrong phi1 and phi2 angle

4 participants